* Analysis of leavers vs remainers
preserve

local bw 150
local donut 0 
local log_y "" // Change name of table below to reflect this variable
local folder `1' // rd_regs

local ef 3.16
local scc_low 50 // central value from IWG 2016
local scc_high 200 // representative value from recent estimates (daniel et al 2019, hansel 2020)

local report_extra_mean 0
local pre_posts 1
local years 2009 2011 
local cases remain switch
local table_tags `2'  // main_supp_levels main_levels  

* update labels, b/c changed in some procedures above
label var dist_eca2009 "ECA Dist, (km)"
label var f_eca09_cons "ECA Fuel, (t)"
label var dist_eca2011 "ECA Dist, (km)"
label var f_eca11_cons "ECA Fuel, (t)"
label var cost_pre_cons "Fuel Cost, No Switch, (USD)"
label var td_inm_pre_cons "Damage, No Switch, (USD)"

local vce "cluster ves_id"


foreach tag of local table_tags {

	foreach year of local years {
		if `year'==2009 {
			local run t_eca09
			local cut eca1
			local y_tag "09"
			local SBcond_pre 1
			*local yvars dist_eca2009 f_eca09_cons dist f_cons cost_comply_cons td_inm_comply_cons cost_eca09_cons td_inm_eca09_cons // tkm_cons
		}
		if `year'==2011 {
			local run t_eca11
			local cut eca2
			local y_tag "11"
			local SBcond_pre 0 
			*local yvars dist_eca2011 f_eca11_cons dist f_cons cost_comply_cons td_inm_comply_cons // cost_eca11_cons td_inm_eca11_cons // tkm_cons
		}

		* getting mean fuel prices for predictions
		* days pre ECA
		local t_opts 30 0
		foreach t of local t_opts {
			sum LS if `run'==-`t'
			local LSmean`t' = r(mean)
			sum IFO if `run'==-`t'
			local IFOmean`t' = r(mean)
		}		
		
		*** Setting Y variables 
		if "`tag'" == "main" {
		if `year' == 2009 {
				local yvars dist_eca`year' avg_speed_eca`year' f_eca`y_tag'_cons cost_eca`y_tag'_cons e_pm_eca`y_tag'_cons td_inm_eca`y_tag'_cons 
			}
			else {
				local yvars dist_eca`year' avg_speed_eca`year' f_eca`y_tag'_cons cost_eca`y_tag'_cons e_pm_eca`y_tag'_cons td_inm_eca`y_tag'_cons  
			}
			local log_y ""
			local report_co2_damage 0 
		}
		
		if "`tag'" == "main_supp" {
			local yvars dist f_cons  hours e_pm_dist_eca`y_tag'_cons
			local log_y ""
			local report_co2_damage 0 
		}	
		
		foreach case of local cases {
			foreach prepost of local pre_posts {
				
				* in/out of sb channel post policy 
				capture drop case_samp
				if "`case'"=="remain" {
					gen case_samp = ( (`cut'==0) & (SBind==`SBcond_pre') ) | ( (`cut'==1) & (SBind==`SBcond_pre') ) 
				}
				
				if "`case'"=="switch"  {
					gen case_samp = ( (`cut'==0) & (SBind==`SBcond_pre') ) | ( (`cut'==1) & (SBind!=`SBcond_pre') ) 
				}
					
				capture drop tmp_*
				* flagging and counting those in sample post policy
				* positive if in sample post cutoff
				gen case_samp_post = ( `cut'==1 & case_samp==1 ) 
				egen tmp_cnt_samp = total(case_samp_post) , by(route_id ves_id)
				drop case_samp_post
				
				* Limit Sample to vessels in sample that we observe pre
					* so keep only vessels we observe leaving post policy
				* or only vessels we observe staying post policy 

				gen tmp_pre = (`run'>=-`bw') & (`run'<-`donut') & tmp_cnt_samp>0
				egen tmp_cnt_pre = total(tmp_pre) , by(route_id ves_id)
				gen tmp_post = (`run'<=`bw') & (`run'>`donut') & tmp_cnt_samp>0 
				egen tmp_cnt_post = total(tmp_post) , by(route_id ves_id)

				gen tmp_pre_post = tmp_cnt_post>0 & tmp_cnt_pre>0 
				tab tmp_pre_post
				
				local prepost_opt "inlist(tmp_pre_post,1,0)" 
				if `prepost'== 1 {
					local prepost_opt "inlist(tmp_pre_post,1)"
				}	
					
				local out_file "`folder'/rd_`tag'_`case'_prepost`prepost'_`year'.tex"
				
				*summ built dwt loa draft power_kw USflag if `prepost_opt' & case_samp==1 & sample_RD_regs==1 & vesseltype_regstr=="Container" & route_name=="LALB_SFBay" ///
				*		& ( `run'<=`bw' & `run'>=-`bw' ) & ( `run'<-`donut' | `run'>`donut' )
				
				local apprep replace
				foreach yvar of local yvars {
				
					* generate sample
					capture drop tmp_samp
					gen tmp_samp = (`prepost_opt') & (case_samp==1) & sample_RD_regs==1 & vesseltype_regstr=="Container" & route_name=="LALB_SFBay" ///
						& ( `run'<=`bw' & `run'>=-`bw' ) & ( `run'<-`donut' | `run'>`donut' ) 
								
					* calculating change in total and specifying formula for spillover
					* setting default			
					if ( "`tag'" == "main_levels" ) | ( "`tag'" == "main" ) {	
						local perc_text ", % of no adjust cf., -9999 "
						local spill_text " , Spillover ratio, -9999 " 
						local report_spill 1
					}
					

					if ( "`tag'" == "main" ) | ( "`tag'" == "main_supp" )  {	
						local perc_text 
						local spill_text " , Spillover ratio, -9999 " 
						local report_spill 1
					}
					
					if ( "`tag'" == "main_supp_levels" ) {	
						local perc_text ", % of no adjust cf., -9999 "
						local spill_text  
						local report_spill 0
					}
					
					if ("`tag'" == "main_corr") {	
						local perc_text 
						local spill_text  
						local report_spill 0
					}
					
					if (`report_spill' == 1) { 
						if "`yvar'" == "dist_eca2009" | "`yvar'" == "dist_eca2011" {
							reghdfe dist `cut' c.`run'#i.`cut' IFO LS if tmp_samp==1 ///
							, absorb(ves_id#route_id) vce(`vce') 	
							local delta_tot = _b[`cut']
							local spill_text ", Spillover ratio, ( 1 - ( `delta_tot' / _b[`cut'] ) ) "
						}
						if "`yvar'" == "f_eca09_cons" | "`yvar'" == "f_eca11_cons" {
							reghdfe f_cons `cut' c.`run'#i.`cut' IFO LS if tmp_samp==1 ///
							, absorb(ves_id#route_id) vce(`vce') 
							local delta_tot = _b[`cut']
							local spill_text ", Spillover ratio , ( 1 - ( `delta_tot' / _b[`cut'] ) ) "
						}
						if ("`yvar'" == "e_pm_comply_cons_ineca09" ) | ("`yvar'" == "e_pm_comply_cons_ineca11") {
							reghdfe e_pm_comply_cons `cut' c.`run'#i.`cut' IFO LS if tmp_samp==1 ///
							, absorb(ves_id#route_id) vce(`vce') 		
							local delta_tot = _b[`cut']
							local spill_text ", Spillover ratio, ( 1 - ( `delta_tot' / _b[`cut'] ) ) "
						}
					}
					

					if "`yvar'" == "cost_eca09_cons" {
						local perc_text ", % of no adjust cf., ( 100 * ( _b[`cut'] / ( `beta_cost_comply_cons' - _b[`cut'] ) ) )"
					} 
					if "`yvar'" == "td_inm_eca09_cons" {
						local perc_text ", % of no adjust cf., ( 100 * ( _b[`cut'] / ( `beta_td_inm_comply_cons' - _b[`cut'] ) ) )"
					}
					if "`yvar'" == "e_pm_eca09_cons" {
						local perc_text ", % of no adjust cf., ( 100 * ( _b[`cut'] / ( `beta_e_pm_comply_cons' - _b[`cut'] ) ) )"
					}
					
					* not reporting percent for 2011
					if `year' == 2011  {
						local perc_text 
					}
					
					if ( "`tag'" == "main" ) | ( "`tag'" == "main_supp" )  {	
						local perc_text 
					}
				
					*  regression
					reghdfe `log_y'`yvar' `cut' c.`run'#i.`cut' IFO LS if tmp_samp==1 ///
						, absorb(ves_id#route_id) vce(`vce') 
					estimates store regresults
					local beta_`yvar' = _b[`cut']	
					capture drop reg_samp
					gen reg_samp = e(sample)	
						
					* predicting outcome at t=0 and t=-30
					margins if reg_samp & `cut'==0 , predict(xb) at( `run'=( -30 ) LS=(`LSmean30') IFO=(`IFOmean30') ) /// 
																  at( `run'=( 0 ) LS=(`LSmean0') IFO=(`IFOmean0') )  post
					*margins, coeflegend

					local mean_t30 = _b[1bn._at]
					local mean_t0 = _b[2._at]	
								
					if `report_extra_mean'==1 {
						local mean_text ", Mean (t=0) , `mean_t0', Mean (t=-30), `mean_t30'"
					}
					else {
						local mean_text ", Mean (t=0) , `mean_t0'"
					}
								
					estimates restore regresults	
						
					* pre ECA mean
					*sum `yvar' if e(sample) & `cut'==0 
					*local mean = r(mean)
					
					if "`log_y'"=="l" {
						* This won't work 
						local chng_text ", Level change , ( `mean' * ( exp( _b[`cut']) - 1  ) )"
												
					}
					if "`log_y'"~="l" {
						local chng_text ", % change , ( ( _b[`cut'] / (`mean_t0' + 0.000000001 ) ) * 100 )"
					}
					
					* getting rid of spillover and percent text if not significant at 5\%
					if ( (2 * ttail(e(df_r), abs(_b[`cut']/_se[`cut']))) > 0.1 ) & "`spill_text'"~="" {
						local spill_text " , Spillover ratio, -9999 " 
						*local perc_text
					}
					
					* get changes in CO2 damages
					local co2_text ""
					if `report_co2_damage'==1 {
						if "`yvar'"=="f_cons" {
							local co2_text ", USD CO2 Damage (SCC=`scc_low' USD/t) , ( _b[`cut'] * `ef' * `scc_low' ) , USD CO2 Damage (SCC=`scc_high' USD/t) , ( _b[`cut'] * `ef' * `scc_high' )"
							disp "`co2_text'"
						}
					}
					
					* calculating difference between dist and resid outcomes in pre-period
					* only adding if there is variable with _d tag
					* and setting perc_text to empty
					noi capture {
						sum `yvar'_d if `cut'==0 & reg_samp					
						local cf_val = r(mean)
					}
					
					if _rc==0 & r(mean)!=0 { 
						local cf_text ", $\Delta$ no behave , `cf_val' , % of no behave , 100 * _b[`cut'] / `cf_val' "						
					}
					else {
						local cf_text 
					}
					disp "`cf_text'"
										
					outreg2 using `out_file', `apprep' lab keep( `cut' ) tex(frag) nonotes nocon /// 
					addstat(Vessels , e(N_clust) `mean_text'  `chng_text' `spill_text' `perc_text' `cf_text' `co2_text' )  
					local apprep append
					
					*tab SBind `cut' if e(sample)	
					
					drop tmp_samp
					
				}
				
				* saving sample indicator for summary stats
				*capture drop samp_SB`SBcond_pre'_prepost`prepost_opt'
				local var "samp_`case'_prepost`prepost'_`year'"
				display "`var'"
				capture drop `var'
				gen `var' = e(sample)
				
			}
			
		}
	}
}

* compare characteristics of switchers
local yvars built dwt loa draft power_kw USflag obs_tot
foreach year of local years {
	capture drop samp_SB_prepost1_`year'
	gen samp_SB_prepost1_`year' = 1 if samp_switch_prepost1_`year'==1
	replace samp_SB_prepost1_`year' = 0 if samp_remain_prepost1_`year'==1
	if `year'==2009 {
		label var samp_SB_prepost1_`year' "Avoider"
	}
	if `year'==2011 {
		label var samp_SB_prepost1_`year' "Returner"
	}
	
	* calculate total obs by vessel
	capture drop obs_tot
	egen obs_tot = count(dist) if samp_SB_prepost1_`year' ~= . , by (ves_id)
	label var obs_tot "Voyages"

	*tabstat `yvars' , by(samp_SB_prepost1_`year')

	local apprep replace
	local out_file "`folder'/rd_leave_remain_balance_`year'.tex"
	foreach yvar of local yvars {
		reg `yvar' samp_SB_prepost1_`year'  , robust 
							
		outreg2 using `out_file', `apprep' lab tex(frag) nonotes 
		local apprep append
		
	}
	
}

restore

/*
* Looking at sample ******************
capture drop temp_samp 
gen temp_samp = (sample_RD_regs==1) & (vesseltype_regstr=="Container") & (route_name=="LALB_SFBay") ///
				& t_eca09>=-`bw' & t_eca09<=`bw'

capture drop obs_pre obs_tot pre_flag 
gen pre_flag = (t_eca09<0)
egen obs_pre = sum(pre_flag) if temp_samp , by (ves_id)	
egen obs_tot = count(dist) if temp_samp , by (ves_id)				

capture drop frequent		
gen frequent = obs_pre>=3		

capture drop voyage_num
sort ves_id temp_samp eca1 date
by ves_id temp_samp eca1 : gen voyage_num = _n if temp_samp & eca1>=1

replace voyage_num = 0 if temp_samp & eca1<1		

lgraph SBind voyage_num if voyage_num<5 , by(dirnw)
tab voyage_num USflag
asdsad	
		
lgraph SBind my if temp_samp , by(USflag)
lgraph SBind my if temp_samp , by(dirnw)
lgraph dist my if temp_samp , by(frequent)
ads

lgraph SBind t_eca09 if temp_samp & obs_tot>4 & dirnw, by(ves_id) ///
legend(off)

* pre ECA patterns by vessel
lgraph SBind t_eca09 if temp_samp & USflag & eca1<1 , by(ves_id) ///
legend(off)
* All US Flagged Vessels use channel pre-ECA

reg SBind USflag if temp_samp & eca1==1
reg SBind dirnw if temp_samp & eca1==1
reg SBind frequent if temp_samp & eca1==1

tab flag_mid SBind if temp_samp & eca1<1 & obs_tot>4 

tab flag_mid SBind if temp_samp & eca1==1 & obs_tot>4 

asdasd
*/
